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ABSTRACT 

We present new high spatial resolution (< Of.'l) 1-5/im adaptive optics images, interferometric 
1.3 mm continuum and 12 CO 2-1 maps, and 350/im, 2.8 and 3.3 mm fluxes measurements of the 
HV Tau system. Our adaptive optics images unambiguously demonstrate that HV Tau AB-C is a 
common proper motion pair. They further reveal an unusually slow orbital motion within the tight 
HV Tau AB pair that suggests a highly eccentric orbit and/or a large deprojected physical separation. 
Scattered light images of the HV Tau C edge-on protoplanetary disk suggest that the anisotropy of 
the dust scattering phase function is almost independent of wavelength from 0.8 to 5/mi, whereas 
the dust opacity decreases significantly over the same range. The images further reveal a marked 
lateral asymmetry in the disk that does not vary over a timescale of 2 years. We further detect a 
radial velocity gradient in the disk in our 12 CO map that lies along the same position angle as the 
elongation of the continuum emission, which is consistent with Keplerian rotation around an 0.5-IMq 
central star, suggesting that it could be the most massive component in the triple system. To obtain 
a global representation of the HV Tau C disk, we search for a model that self-consistently reproduces 
observations of the disk from the visible regime up to millimeter wavelengths. We use a powerful 
radiative transfer model to compute synthetic disk observations and use a Bayesian inference method 
to extract constraints on the disk properties. Each individual image, as well as the spectral energy 
distribution, of HV Tau C can be well reproduced by our models with fully mixed dust provided 
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grain growth has already produced larger-than-interstellar dust grains. However, no single model can 
satisfactorily simultaneously account for all observations. We suggest that future attempts to model 
this source include more complex dust properties and possibly vertical stratification. While both grain 
growth and stratification have already been suggested in many disks, only a panchromatic analysis, 
such as presented here, can provide a complete picture of the structure of a disk, a necessary step 
towards quantitatively testing the predictions of numerical models of disk evolution. 
Subject headings: planetary systems: protoplanetary disks - stars: pre-main sequence - stars: individ- 
ual (HV Tau) 



1. INTRODUCTION 

Circumstellar disks are an ubiquitous outcome of the 
stellar formation process and they are believed to be the 
birth place of planetary systems. The growth of dust par- 
ticles towards planetesimal sizes along with their vertical 
settling due to gas drag are processes that are believed 
to be the first steps towards planet formation. Hydrody- 
namical models have shown that these processes can be 
efficient early in the disk evolution (e.g., Weidenschilling 
1997; Dullemond & Dominik 2004; Barriere-Fouchet ct 
al. 2005). To test these models, it is necessary to ob- 
tain an observation-based description of the structure 
and dust content of protoplanetary disks as a function 
of the age of the system and other relevant parameters 
(e.g., stellar mass). 

The dust component of protoplanetary disks has long 
been studied via its thermal emission from near-infrared 
to millimeter wavelengths which is frequently associated 
with low-mass pre-main sequence T Tauri stars (Kenyon 
& Hartmann 1987; Bertout et al. 1988; Strom et al. 1989; 
Beckwith et al. 1990). Both grain growth and dust set- 
tling can alter the overall shape of the spectral energy dis- 
tribution (SED) of a T Tauri star (D'Alessio et al. 2001; 
Dullemond & Dominik 2004; D'Alessio et al. 2006). In- 
deed, several studies that analyzed (elements of) the SED 
of young stars have concluded that both grain growth and 
settling is occuring in protoplanetary disks (e.g., Beck- 
with & Sargent 1991; Mannings & Emerson 1994; Furlan 
et al. 2006; Kessler-Silacci et al. 2006; Rodmann et al. 
2006; Natta et al. 2007, and references therein). Un- 
fortunately, such studies suffer from the absence of spa- 
tial information inherent to photometric measurements 
and the high optical thickness of disks in the near- to 
mid-infrared regime. As a result, comparing an object's 
SED to radiative transfer models leaves many ambigui- 
ties (Chiang et al. 2001). For instance, the inferred total 
dust mass and the maximum size of the dust grains are 
inversely correlated because of the dependency of dust 
opacity on grain size. In addition, most of these stud- 
ies, which focus on a single type of observations (e.g., 
millimeter fluxes, silicate emission feature), only probe 
a limited region of the disk and a small fraction of the 
entire grain size distribution. 

To solve for the ambiguities inherent to SED stud- 
ies, it is critical to obtain spatially-resolved observa- 
tions. Such observations include thermal emission map- 
ping with (sub)millimeter interferometers (e.g., Keene 
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& Masson 1990; Simon ct al. 1992; Lay et al. 1994; 
Dutrey et al. 1996; Andrews & Williams 2007) and scat- 
tered light imaging with optical and near-infrared high- 
resolution instruments (e.g., Burrows et al. 1996; Rod- 
dier et al. 1996; Stapelfeldt et al. 1998). Disks are gen- 
erally optically thin at long wavevengths, so the former 
type of observations can probe the entire disk structure. 
Furthermore, they are very sensitive to the presence of 
millimeter-sized particles. On the other hand, scattered 
light images, which only probe dust grains at the disk 
surface, are very sensitive to the size distribution of mi- 
cronic grains especially when images at multiple wave- 
lengths are analyzed simultaneously (Watson et al. 2007, 
and references therein). Both types of observations have 
already yielded important pieces of evidence support- 
ing both grain growth and dust settling in disks (e.g., 
Duchcnc ct al. 2003, 2004; Watson & Stapelfeldt 2004). 

While grain growth and settling appear to occur in 
protoplanetary disks, detailed quantitative tests of hydro- 
dynamical models can only be achieved with a detailed 
view of the entire structure of a disk. This can only be 
obtained via a multi-technique, panchromatic approach. 
Unfortunately, observational and computational limita- 
tions have so far limited the number of objects for which 
such an analysis could be conducted to a handful. The 
most notable examples are the studies of the "Butterfly 
Star" (Wolf et al. 2003), IMLup (Pinte et al. 2008) and 
IRAS 04158+2805 (Glauser et al. 2008). In the former 
two cases, these studies have unambiguously shown that 
the dust population is stratified, possibly indicating that 
dust settling is already occurring. Increasing the number 
of disks studied in such detail is necessary to disentangle 
individual peculiarities from genuine trends associated 
with disk evolution. 

HV Tau is a triple system located in the Taurus star- 
forming region. It consists of a 550 AU wide pair whose 
optically brightest component is itself a tight (10 AU) 
visual binary (Simon et al. 1996). Spectroscopic and 
photometric measurements revealed that this subsystem 
does not currently experience accretion nor does it show 
infrared excess. They further establish an age of about 
2Myr for the system (White & Ghez 2001). The third 
component of the system, HV Tau C, is much fainter yet 
bluer than HV Tau AB. While Magazzu & Martin (1994) 
first thought that this source was an Herbig-Haro object, 
Woitas & Leinert (1998) later proposed that HV Tau C is 
a normal M0 T Tauri star surrounded by an opaque edge- 
on disk similar to that found in HH 30 by Burrows et al. 
(1996). Subsequent high resolution imaging confirmed 
this hypothesis (Monin & Bouvier 2000). Stapelfeldt et 
al. (2003) produced the first model of high resolution 0.8 
and 2.2^m scattered light images of HV Tau C, finding 
that dust properties similar to those of interstellar dust 
grains can account for these images. Their 0.8^m image 
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also revealed the presence of a roughly spherical enve- 
lope, producing a symmetric halo that is more extended 
than the disk itself, that is likely the remnant of the core 
from which the system formed. 

As a consequence of their particular viewing geometry, 
edge-on protoplanetary disks offer a unique opportunity 
to determine their geometry and dust content. As such, 
they may be the best candidates to study vertical strat- 
ification in protoplanetary disks. They also are compar- 
atively easy targets for high-angular resolution instru- 
ments since the contrast requirement is strongly relaxed. 
On the other hand, they are challenging from the mod- 
eling point of view because of the difficulty of accurately 
solving for radiative transfer including anisotropic scat- 
tering in high optical depth regions. Previous modeling 
efforts of edge-on protoplanetary disks have therefore fo- 
cused on interpreting one type of observation at a time 
(Burrows et al. 1996; Stapelfeldt et al. 1998, 2003; Cotera 
et al. 2001; Wood et al. 2002; Watson & Stapelfeldt 
2004). While these studies proved highly valuable to 
constrain some of the disk properties, no self-consistent 
model was used to model all data at once, leaving unex- 
plained contradictions. 

Our objective in this work is to perform a global anal- 
ysis of the HV Tau C disk, combining scattered light 
images, millimeter interferometric data and the overall 
SED into a single fit. We present a series of new observa- 
tions of the system in Section 2 and discuss our empirical 
results in Section 3. In Section 4, we present radiative 
transfer models of the HV Tau C disk and discuss their 
implications in Section 5. Section 6 summarizes our re- 
sults. 

2. OBSERVATIONS 
2.1. Adaptive optics near-infrared imaging 
2.1.1. l-2^m imaging 

On 2002 November 24, we observed HV Tau using the 
NAOS adaptive optics system and the CONICA instru- 
ment (Lenzen et al. 2003; Rousset et al. 2003) installed 
on the Yepun 8.2 m-Unit Telescope at ESO's Very Large 
Telescope, as part of the NAOS Guaranteed Time Ob- 
serving program. We used the 0'.'0133 pixel scale for the 
J (A = 1.27^m, AA = 0.25/zm) and H (A = 1.66^m, 
AA = 0.33/im) images and the 0'.'0270 plate scale for the 
K s (A = 2.18^m, AA = 0.35/im) images. HV Tau AB, a 
my = 14.5 source, was used as adaptive optics guide star 
with the visible wavefront sensor. From narrow-band ob- 
servations of single stars throughout the night, the mea- 
sured FWHM of point-like sources is about O'.'IO, 0"07 
and 0'.'08 at J, H and K s band, respectively. These im- 
ages, similarly, have Strehl ratios of approximately 5%, 
25% and 40%. 

We obtained deep images in which HV Tau AB is sat- 
urated, as well as shallow images with an additional neu- 
tral density filter and shorter exposures to record unsat- 
urated images of the primary. Total integration times for 
the long exposures of 750s, 480s and 165.5s were recorded 
at J, H and K s , respectively, split in 8 to 20 dithered 
independent images. Total integrations times for the 
shallow images were 100s, 50s and 45s, respectively at 
J, H and K s . For each sequence of images, a sky was 
estimated by medianing all images, which was then sub- 
tracted from each image prior to cosmetic cleaning, which 



included bad pixels and cosmic ray correction, and flat- 
fielding. All resulting images were then shift-and-added 
to produce final images. 

On 2002 November 26, we re-observed HV Tau with 
NAOS and CONICA using the 0'.'0270 pixel scale, this 
time using the l'.'4-diameter coronagraph mask to block 
out the starlight from HV Tau AB and record deeper 
exposures on HV Tau C. Three 120s K s band images 
were recorded. A sky was subtracted off each image be- 
fore they were cosmetically cleaned and flat-fielded. Fi- 
nally, they were averaged to yield the final coronagraphic 
image. No point source is detected in the 30" field-of- 
view, preventing us to estimate the achieved FWHM and 
Strehl ratios. 

2.1.2. 3-5[im imaging 

On 2002 December 13, we observed HV Tau using 
the NIRC2 camera (P.I. K. Matthews) installed behind 
the Keck II adaptive optics system (Wizinowich et al. 
2000) to record images using the L' (Ao = 3.78/im, 
AA = 0.70^m) and M s (A = 4.67yum, AA = 0.24^m) 
filters and the 0'.'00996 plate scale (Ghez et al. 2008). 
HV Tau AB was used as a Natural Guide Star (NGS) for 
the adaptive optics system. To reduce thermal back- 
ground and provide a more symmetric and smoother 
PSF, we used the "inscribed circle" pupil stop, result- 
ing in an effective primary mirror diameter of 9m. At 
L', four images, consisting of 300 coadded 0.181s indi- 
vidual integrations each (for a total of 217.2s on source 
integration), were acquired with the sources dithered on- 
chip between images. A fifth image was acquired with 
all stars moved out of the detector field-of-view (10") 
to estimate the sky level. At M s , 500 individual 0.150s 
integrations were coadded using a reduced 6"-wide field- 
of-view, necessary to avoid background saturation on the 
intense thermal background. A total of 21 such images 
were obtained with the stars moved about in the avail- 
able field-of-view, representing a total integration time of 
1575s. The corresponding sky was estimated by median- 
combining all frames. 

On 2004 November 03, we re-observed HV Tau with 
NIRC2, this time using the newly available Laser Guide 
Star (LGS) module on the Keck II adaptive optics system 
(van Dam et al. 2006; Wizinowich et al. 2006). The laser 
system was run at 6 Watts output power, which produced 
a guide star with an equivalent magnitude of my = 9.9, 
and HV Tau AB was used as the tip/tilt correction point 
source. To take full advantage of the better image qual- 
ity, the "largehex" pupil mask, which does not block any 
section of the primary mirror, was used throughout the 
observations. HV Tau was imaged in the L' filter with 
an 0.2s integration coadded 100 times. This observing 
cycle was repeated 27 times, with the system alterna- 
tively located in opposite detector quadrants, providing 
a total integration time of 540 sees. For each image, 
the subsequent one, with the star located in the opposite 
quadrant, was used as a sky frame. 

All datasets were reduced using a similar strategy 
to that used in processing the NACO images (see also 
Duchcnc et al. 2004). First of all, the sky thermal emis- 
sion was subtracted. These subtracted images were then 
flat fielded and had any bad pixels interpolated over. In 
the case of the LGS dataset, a residual sky level in the 
frames was measured by taking the median value of the 
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quadrants which do not have a star in the field. This 
value (typically corresponding to 0.05% of the initial sky 
value) was then subtracted off the cleaned images. All 
cleaned images in a given dataset were then median- 
averaged. We measured FWHMs of 0'.'09 and O'.'ll for 
point- like sources observed before/after HV Tau at L 1 
and M s for the NGS dataset and 0"08 at L' for the LGS 
dataset. We estimated the corresponding Strehl ratios to 
be about 70% at L' (both in NGS and LGS mode) and 
80% at M s , using the "Strehl meter" tool developed by 
the Keck observatory 2 . 

Given the high thermal background at L' and M s , the 
resulting signal-to-noise on HV Tau C is limited, about 
15, 20 and 5 in the peak pixel at V (NGS, LGS) and M s , 
respectively. To enhance detection without altering the 
intrinsic spatial resolution of the images, we smoothed 
them by using a running 2-pixel-radius (~0'.'02) median- 
filtering circular mask, and rebinned them by a factor 
of 3 in both directions, resulting in an approximate final 
sampling of 0703/pixel that still oversamples the resolu- 
tion of the datasets. 

2.2. Sub-millimeter and millimeter imaging 
2.2.1. 1.3 and 2.8 mm interferometric imaging 

On 2005 February 26, we observed HV Tau with 
IRAM's 6-antcnna Plateau de Bure Interferometer 
(PdBI, Guilloteau et al. 1992) in the 6Bp configuration, 
with baselines ranging from 71m to 331m. Simulta- 
neous 110 GHz (2.76 mm) and 230.5 GHz (1.31mm) ob- 
servations of HV Tau were recorded in double-sideband 
mode with a bandpass of 640 MHz at each frequency. 
HV Tau was observed alternatively with phase calibra- 
tors 0415+379 and 0528+134 through an entire llh- 
transit, resulting in beam sizes of l'.'lxO'/9 (along po- 
sition angle 38°) and 2'.'lxl'.'6 (along position angle 62°) 
at 1.3 and 2.8 mm, respectively. The average weather 
conditions resulted in r.m.s. phase noises on the order 
of 15-40° at 2.8 mm and 8-15° at 1.3 mm (using self- 
calibration from the 2.8mm data), equivalent to an at- 
mospheric "seeing" on the order of Of! 3 and C/6 at 1.3 
and 2.8 mm, respectively. The absolute pointing uncer- 
tainty is on the order of 0"l-0"2. The quasars NRAO 150 
and 3C 273 were used as absolute flux calibrators, result- 
ing in a 10% uncertainty on all quoted fluxes. The data 
were reduced using the GILDAS package and selecting 
individual baseline visibilities for which the phase noise 
was less than 40° and flux variations were less than 20% 
based on calibrators measurements. Simultaneous obser- 
vations of the 12 CO 2-1 transition (230.538 GHz rest fre- 
quency) with a 20 MHz bandpass were obtained to probe 
the gaseous component of the disk. Data reduction fol- 
lowed the same method as the continuum data, providing 
a 3-dimension reconstructed datacube with a 0.1 km.s 
spectral resolution. 

We obtained follow-up observations of HV Tau with 
the 15-antenna CARMA array in its C configuration, 
with baselines ranging from 24 to 300m. On 2008 May 
6 and 9, we tuned the receivers to a central frequency 
of 110 GHz (2.76mm), while on 2008 May 29 we tuned 
them to 90 GHz (3.33 mm). The total continuum band- 
width is 2.8 GHz, split in 6 separate bands. Observing 

2 http://www2.keck.hawaii.edu/optics/aochar/Strehl_meter2.htm 



conditions were average and parts of the observations 
had to be flagged out because of poor phase coherence. 
Overall, the useful integration times on HV Tau were 
3h45 at 110 GHz and 2h30 at 90 GHz. Observations of 
HV Tau were interleaved with pointings at 3C111 and 
J 0530+135 which served as phase calibrators; flux cal- 
ibration was performed by observing a planet (Uranus, 
Neptune) at the beginning of each track. The systematic 
uncertainty in CARMA's absolute flux scale is ~20% (W. 
Kwon, priv. comm.). The data were reduced using the 
MIRIAD software package. The final 110 GHz map cor- 
responds to the combination of both observing periods, 
and is characterized by a beam size of l'/8xl'/4 (along 
position angle 116°); the 90 GHz map has a beam size of 
3'.'lxT/7 (along position angle 124°). 

2.2.2. 350[im single dish mapping 

On 2008 January 28, we observed HV Tau at 350^m 
with the 32x12 SHARC-II bolometer array (Dowell et 
al. 2003) installed at the Caltech Submillimeter Observa- 
tory. The array was scanned at an 8'.'2/s rate over a 60" x 
field centered on HV Tau. This scanning was repeated 
until a total integration of 600 s was achieved. Ceres and 
HLTau, two secondary flux calibrators, were observed 
immediately after with total integration times of 180 s 
and 300 s, respectively. The fluxes of both these sources 
are known to within 10%; their measured fluxes agreed 
to within this uncertainty. Conditions were excellent for 
submillimeter observations (0.013 < T225GHZ < 0.036). 
All data were reduced using the CSO-developed CRUSH 
software, using a 6" smoothing for HV Tau (a 4" smooth- 
ing was used for the flux calibrators), to produce a final 
image that has an effective resolution of approximately 
10". 

3. OBSERVATIONAL RESULTS 
3.1. HV Tau as a triple system 

At all wavelengths, our adaptive optics images clearly 
reveal the typical morphology of HV Tau C as an opaque 
edge-on disk, namely two parallel, horizontal nebulae 
separated by a dark lane (Figure 1). They further show 
HV Tau AB to be systematically extended along the 
north-western direction, although the pair is barely re- 
solved due to its tight separation (see Figure 2). Here we 
discuss the relative astrometry of all 3 components. 

Focusing first on the wide pair HV Tau C-HV Tau AB, 
we estimate its separation and position angle based on 
the location of the centroid of both components. De- 
spite saturation of HV Tau AB in some images, all im- 
ages yield consistent estimates (see Table 1), with an 
average separation of 4'.'04±0"02 and a position angle 
of 44?6±0?6. There is marginal (« 3a) evidence for a 
variation of the binary position angle as a function of 
wavelength, though this could be due to underestimated 
uncertainties on the absolute orientation of the detectors. 
We do not find any significant trend as a function time 
at our level of precision. Based on the proper motion of 
HV Tau AB measured by Ducourant et al. (2005), the 
projected separation and position angle of the pair would 
have changed in the almost 5 years between the obser- 
vations of Stapclfeldt et al. (2003) and ours by a total 
of 0"12 and 2?2, respectively, if HV Tau C had not been 
comoving with HV Tau AB. This can be rejected at the 
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« 5cr level in our datasets, and is further confirmed by the 
relative astrometry obtained at earlier epochs by Simon 
et al. (1992) and Woitas & Leinert (1998). We therefore 
conclude that HV Tau is a common proper motion pair 
with a projected separation of 565 AU, making it a bona 
fide triple system. 

To obtain reliable relative photometry and astrometry 
for the tight HV Tau AB pair, we used PSF-fitting. We 
choose to analyze the short H- and if s -band VLT expo- 
sures and both L'-band Keck images, as they offer the 
best image quality and the most favorable ratio between 
binary separation and achieved resolution. For the VLT 
-ff-band and Keck L'-band images, we searched for ad- 
equate PSFs among images of single stars obtained on 
the same night as our HV Tau images. For the VLT 
.Ks-band, we used a set of single stars from observations 
taken with the same set-up presented in Duchene et al. 
(2007). PSF-fitting was performed using IRAF's daophot 
package. The resulting relative astrometry and photom- 
etry is listed in Table 1. The VLT and Keck NGS data 
have been taken within a few weeks of each other, and 
no measurable orbital motion is expected over such a 
short timescale. Averaging these three datasets, we find 
a separation of 0'.'0595±0'.'0025 and a position angle of 
312?5±1?8. The Keck LGS data, obtained two years 
later, yields a relative astrometry that is marginally con- 
sistent (2.5<j) with this estimate. The projected sepa- 
ration in our images appears to be 4-5ct smaller than 
earlier measurements (Simon et al. 1996; Monin & Bou- 
vier 2000), but there is no significant change in position 
angle as a function of time (see Figure 2). This rela- 
tive displacement is most likely due to orbital motion. 
The observed plane-of-the-sky velocity of HV Tau AB, 
«1.5km.s _1 , is roughly one order of magnitude too low 
considering the 10 AU projected separation of the bi- 
nary compared to other TTauri binaries (e.g., Ghez et 
al. 1995). The number of resolved measurements and 
the total amplitude in orbital motion are insufficient to 
attempt an orbital fit for HV Tau AB at this point. 
Nonetheless, the unexpectedly low measured orbital ve- 
locity suggests a highly elliptical orbit observed around 
apoastron passage, or a large out-of-the-plane separation 
which, combined with the nearly radial observed motion, 
would in turn imply that the orbital plane is almost per- 
pendicular to the plane of the sky. Monitoring of the 
system in the next few years will help disentangle these 
two possibilities. 

Our PSF fitting also yields relative photometry for the 
HV Tau AB binary. Considering the near-simultaneous 
H, K s and NGS L' images, we find evidence that 
HV Tau B is somewhat redder than HV Tau A in the 
near-infrared (see Table 1). In the framework in which 
none of the components possesses circumstellar mate- 
rial, this is an indication that HV Tau B is cooler than 
HV Tau A, consistent with it being fainter. Surprisingly, 
Simon et al. (1996) found a flux ratio in the visible that 
is closer to unity than our near-infrared measurements. 
Temporal variability of either component could be an ex- 
planation, although we note that both our L' flux ratios 
are consistent with one another despite being taken two 
years apart. Further monitoring is required to better 
understand the intrinsic colors of both components. 

3.2. The HV Tau C circumstellar disk 



3.2.1. Scattered light images 

The near-infrared images presented here are of higher 
spatial resolution than those presented in Monin & Bou- 
vier (2000) and Stapelfeldt et al. (2003), which both had 
a 0'.'13 resolution, and are comparable to the HKL 1 Sub- 
aru images of Tcrada et al. (2007). Our observations ex- 
tend the wavelength coverage of HV Tau C to M s for 
the first time and, combined with the HST F814W im- 
age of Stapelfeldt et al. (2003), offer an almost uniform 
spatial resolution (O'.'07-O'.'ll, or 10-15 AU) view of the 
disk over the entire 0.8-4.7/xm range. This is critical to 
conduct unbiased studies of the wavelength dependence 
of the disk images. 

To quantify the basic morphological properties of 
HV Tau C, we adopt the following method. At each 
wavelength, we first estimate the position angle of the 
dark lane, PA^isk, as the average of the position angles of 
each nebulae, as determined from fits of elliptical Gaus- 
sian intensity profiles. Averaging all estimates, we find 
a position angle of 108?3±0?4, which we take as the ori- 
entation of disk midplane. The dark lane width, d nc b, 
is measured as the projection of the vector joining the 
light centroid of the two nebulae on the disk minor axis. 
Peak-to-peak flux ratios (-Fi? p0 ak) arc readily estimated, 
whereas integrated flux ratios (FRi nt ) are obtained by 
summing the flux within areas that encompass all pixels 
whose surface brightness is at least 5% (20% for the M s 
image) of the peak surface brightness in the image. Fi- 
nally, we measure the total extent of the disk along its 
major-axis, w^ a , defined by the horizontal extent of the 
contour at the 5tr noise level. We note that the spherical 
halo identified by Stapelfeldt et al. (2003) dominates this 
measurement in the F814W image. All of these quanti- 
ties are given in Table 2. 

Several wavelength-dependent features can be noted in 
our images of HV Tau C. As Figure 3 illustrates, c? ne b de- 
creases by about 30% from 0.8 to 4.7^m. We further find 
that both Fi? pca k and FR- lnii smoothly decline towards 
longer wavelengths without ever reaching unity; rather, 
both flux ratios actually plateau longwards of 2.2/um. Fi- 
nally, the intensity profile of the bright nebula along the 
disk major axis is remarkably invariant from 0.8 to 4.7/im 
(see Figure 4). The measured FWHM is about 0'.'35 at all 
wavelengths, i.e., a factor of 3-4 larger than the resolu- 
tion of our images. This achromatic behavior is therefore 
not influenced by observational limitations but is rather 
intrinsic to the disk. 

One main feature of the HV Tau C disk that was 
pointed out by Stapelfeldt et al. (2003) is the lateral 
asymmetry: the northeastern (bright) nebulae extends 
further on one side while the southwestern (faint) nebula 
extends further on the other side. This can be read- 
ily seen in Figure 3. In addition to this low-intensity 
asymmetry, the location of the nebulae ccntroids is also 
non symmetric, suggesting a global asymmetry: the line 
joining the two centroids is misaligned by about 10° with 
respect to the disk symmetry axis as determined by the 
orientation of the dark lane. We do not find any sig- 
nificant variation of the asymmetry with wavelength nor 
time within our uncertainties. 

Stapelfeldt et al. (2003) also discussed narrow struc- 
tures, which they dubbed "rays", extending along the 
disk major axis beyond the main/bright nebula which 
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they associate to the disk, with total extent up to l'/4 in 
their optical images. We can track scattered light up sim- 
ilar distances at both K s and L' in our images, a regime 
in which the halo is undetected. Although the exact na- 
ture of these features remains uncertain, this suggests 
that these "rays" trace the surface of the disk instead of 
being a mere shadow of the disk on the spherical halo. 
In that case, the actual disk size may be larger than the 
50 AU assumed in the past. 

3.2.2. Thermal emission regime 

Both our PdBI and CARMA observations have a suf- 
ficient spatial resolution to disentangle the wide HV Tau 
pair. In all maps, a single source is detected, at the lo- 
cation of HV Tau C. HV Tau AB appears to have no 
significant emission at millimeter wavelengths with Sc- 
upper limits in our PdBI data of 6mJy and 1.5 mJy at 
1.3 and 2.8mm, respectively. 

At 1.3 mm, HV Tau C is spatially resolved in our PdBI 
data, with a sharp decrease in correlated flux at the 
longest baselines along the position angle of the disk as 
defined from the scattered light images (see Figure 5). 
We therefore fit an elongated Gaussian profile to the 
visibilities, assuming that its minor axis is unresolved 
which is consistent with the data. This yields a total 
flux density of 49.5±1.8mJy, a FWHM of l'/17±0'/08 
along the major axis, at PA 111±3° East of North, in 
excellent agreement with the position angle derived from 
the scattered light images. This analytical fit is shown 
in Figure 5; deviations between observations and our 
simple fit result from the fact that the distribution of 
surface brightness in the disk is not a perfect Gaussian. 
The measured 1.3 mm flux agrees well with the previous 
single-dish measurement by Beckwith et al. (1990). At 
2.8 mm, the source is unresolved with the PdBI, with a 
flux density of 7.1±0.5mJy. The CARMA observations 
do not resolve the disk either, and we extracted point- 
source fluxes of 8.0±2.1mJy at 2.8 mm and 3.8±0.9mJy 
at 3.3 mm. Flux calibration uncertainties of 10% and 
20% must be added for the PdBI and CARMA observa- 
tions, respectively. 

A single point-like source is detected in our 350^m 
CSO map of HV Tau with a flux of 0.370±0.030 Jy; a 
10% uncertainty for flux calibration must be quadrati- 
cally added. The resolution of this dataset, roughly 10", 
does not allow us to spatially resolve the wide pair. How- 
ever, taking advantage of the strong detection of the sys- 
tem (« 20er) , we use the nearly simultaneous observation 
of the bright source HL Tau with the same set-up to de- 
termine an accurate astrometric positioning (Mundy et 
al. 1996). The 350£im source is found to be located at 
(04h38m35.50s, +26°10'40"6, J2000) with an uncertainty 
of about 1" in both directions. This is in excellent agree- 
ment both with the position of the lone millimeter source 
in our PdBI and CARMA maps and with the location of 
HV Tau C in optical and near- infrared images. We there- 
fore conclude that HV Tau C is the dominant source of 
emission at 350^im and assign all of the observed flux to 
that component. 

3.2.3. Gas emission 

As shown in Figure 6, 12 CO 2-1 line emission is de- 
tected at the position of HV Tau C in our PdBI observa- 
tions. The integrated spectrum shows two distinct peaks, 



with a near-zero trough separating them. The trough 
in the middle of the line is likely due to contamination 
by large-scale emission from the surrounding molecular 
cloud that is filtered out by the interferometer. Mizuno 
et al. (1995) detected 13 CO emission in the 5.5-7km s _1 
from the molecular cloud at the location of HV Tau, con- 
sistent with this interpretation. We further find that the 
blue- and red-shifted parts of the line emission are spa- 
tially distinct and symmetric about the continuum peak. 
The two line emission peaks are separated by about l'/5, 
or 200 AU, spatially and 4km.s _1 in velocity. Finally, 
we note that the CO emission appears to extend beyond 
both the millimeter continuum and the scattered light 
images, suggesting that there is more to the system than 
meets the eye in continuum datasets. 

3.2.4. The SED of HV Tau C 

To draw a complete picture of the HV Tau C disk, 
we compiled its SED by combining published fluxes with 
our new measurements. From the same references (see 
below), we also built the SED of HV Tau AB (Figure 7), 
which is well reproduced from the optical to the mid- 
infrared by a 3600 K, logg = 4.0, NextGen model from 
Baraffe et al. (1998) assuming Ay — 1.75 mag, in agree- 
ment with previous spectrophotometric estimates of the 
stellar properties of that component (Kenyon & Hart- 
mann 1995; White & Ghez 2001). 

Constructing the SED of HV Tau C is a delicate mat- 
ter because of the large scatter in its published optical 
and near-infrared photometry (see also Monin & Bou- 
vier 2000; Stapelfeldt et al. 2003). Figure 7 includes all 
published data and illustrates this variability. To con- 
struct a "representative" SED for HV Tau C, we elected 
to aim for smoothness. For instance, we adopt the most 
recent KLN photometry from McCabe et al. (2006), 
which matches well with the Spitzer/IRAC data from 
Hartmann et al. (2005). We discard the N band mea- 
surements from Woitas & Leinert (1998) as they appear 
to overestimate the flux of both components by about 
1 mag. Based on their agreement at K band with the Mc- 
Cabe et al. photometry, we then adopted the JH fluxes 
from Woitas & Leinert (1998). Finally, we adopted the 
visible photometry from Magazzu & Martin (1994) which 
offers a smoother extension of the SED to short wave- 
lengths than that of Stapelfeldt et al. (2003), although 
it was obtained almost a decade earlier. Apart from the 
K band measurement by Simon et al. (1992), our selec- 
tion is equivalent to consistently selecting the brightest 
measurement available at every wavelength. At longer 
wavelengths, we used the MIPS 24/xm and 70/zm detec- 
tions and a 170/xm upper limit from the Spitzer Taurus 
Legacy Survey (Rebull et al., submitted). At 24^m the 
astrometric accuracy is good enough to identify the de- 
tected source with HV Tau C with little or no contribu- 
tion from HV Tau AB. In the longest wavelength regime, 
we adopt fluxes measured from our data and from An- 
drews & Williams (2005). The resulting SED, the most 
complete to date for an edge-on disk system, is shown as 
filled diamonds in Figure 7 and presented in Table 3. 

The SED of HV Tau C can be characterized as a "dou- 
ble hump" (see Figure 7) , which is typical of other edge- 
on disks (e.g., Strom & Strom 1994; Stapelfeldt & Mon- 
eti 1999; Wood et al. 2002). The short wavelength hump 
is dominated by photons emitted from the central star 
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and the innermost region of the disk and scattered to 
the observer at the outer radius of the disk, whereas the 
long wavelength part represents the disk thermal emis- 
sion propagated through the mostly optically thin disk 
to the observer. The wavelength of the turnover between 
the two humps is driven by the total column density of 
dust along the line of sight to the central star. We find 
this turnover to be between 8 and 11.8/im although a 
more precise estimate would require analysis of a mid- 
infrared spectrum of the source. At longer wavelengths, 
Andrews & Williams (2005) noted that HV Tau C has 
an abnormally flat 850/im-1.3mm slope. With our new 
observations, it appears that the 850/im flux is indeed 
too low compared to extrapolation from longer wave- 
length fluxes. Limiting our analysis to the 1.3-3.3 mm 
regime, wc find a spectral index of a mm = 2.5±0.2, where 
F v oc v a ™ m . This value is on the low end of the distribu- 
tion observed for other protoplanetary disks (Bcckwith 
& Sargent 1991; Mannings & Emerson 1994; Natta et al. 
2007). 

3.2.5. Qualitative interpretation 

Here wc interpret some of the observational results 
listed above to make qualitative inferences about the 
HV Tau C system. We will then revisit in a quanti- 
tative way these conclusions on the basis of the radiative 
transfer modeling presented in Section 4. 

The >lmag variability of HV Tau C in the near- 
infrared is larger than has been recorded for other 
TTauri stars in the past (Eiroa ct al. 2002; Alves de 
Oliveira & Casali 2008). It is therefore unlikely that 
this variability results from intrinsic fluctuations in the 
emission of the star and inner disk. A plausible alterna- 
tive scenario could invoke a self-shadowing pattern in the 
disk that moves as a result of differential Keplerian rota- 
tion. In this framework, the structure responsible for the 
lateral asymmetry in the images would be located a few 
AU away from the central star to account for the absence 
of variability over the course of 2 years. This could be 
linked to the observed asymmetry in the scattered light 
images. Such a phenomenon has been documented in the 
HH30 edge-on disk (Watson & Stapelfeldt 2007). One 
possibility is that the disk is warped as a consequence of 
the gravitational forces exerted by HV Tau AB. 

The wavelength dependent features of the scattered 
light images inform us about the absorption and scat- 
tering properties of the dust grains in the disk. For in- 
stance, the decreasing d ne b towards longer wavelengths is 
a consequence of the declining dust opacity (Stapelfeldt 
et al. 2003; Watson & Stapelfeldt 2004). The fact that 
the intensity profile along the bright nebula does not 
change with wavelength is suggestive of a dust popu- 
lation whose phase function is equally anisotropic at all 
wavelengths. A dust population similar to that of the in- 
terstellar medium, whose scattering asymmetry param- 
eters drops rapidly beyond l^m (e.g., Weingartncr & 
Drainc 2001), seems inconsistent with our observations 
of HV Tau C. We note, however, that the flux ratio be- 
tween the nebulae is not constant across the 0.8-4.7^im 
range, even though it is also believed to be driven by the 
scattering asymmetry. 

The slope of the SED of an edge-on disk between the 
near- and mid-infrared is a function of the disk geom- 
etry and of the grain size distribution and composition 



in the disk scattering surface. In their analysis of the 
"Flying Saucer" edge-on disk, Pontoppidan et al. (2007) 
concluded that the shallow 2-10^m slope in that system 
implies the presence of a "significant amount of 5-10/zm 
grains" in that disk. Such grains are necessary to produce 
a high albedo in the mid-infrared. While it is tempting 
to apply this argument to HV Tau C, uncertainties about 
the disk geometry prevents from reaching a firm conclu- 
sion just yet. 

If the disk were optically thin to its own emission in 
the millimeter regime, the observed spectral index could 
be readily converted into an opacity power law index of 
0.5 ± 0.2, characteristic of evolved dust grain popula- 
tions that extend up to a few millimeter in size (Man- 
nings & Emerson 1994; Natta et al. 2000; Andrews & 
Williams 2005). However, because the disk is very com- 
pact and viewed almost exactly edge-on, the optically 
thin assumption may be incorrect, so that grain growth 
cannot be claimed on this sole basis in this particular 
system. If the disk were indeed optically thick, one 
would expect the millimeter spectral index to steepen 
at longer wavelengths. Wc find a 1.3-2. 7mm = 2.5 ± 0.2 
and 

^2.7— 3.3mm — 

3.4 ± 0.6, a difference that is not sta- 
tistically significant. Flux measurements at even longer 
wavelengths are ultimately needed to conclude on this 
possibility. 

In any case, there are several independent hints of the 
presence of grains at least a few microns in size, and 
maybe up to millimeter sizes, in the HV Tau C disk. 
However, none of these inferences are robust enough as 
they depend on some assumptions about the disk ge- 
ometry, among other factors. The goal of the modeling 
presented in the next section is to simultaneously fit for 
the disk structure and the dust properties so as to solve 
for these ambiguities. 

In the absence of a gas tracer whose emission from 
the cloud is negligible, we do not perform a complete 
Keplerian rotation analysis of our 12 CO observations. 
Qualitatively, the measured amplitude of the velocity and 
positional offsets are consistent with an 0.5-1Mq cen- 
tral star (see Figure 6). Observations in other molecular 
tracers are necessary to go beyond this simple analysis, 
however. We also note that the larger outer disk radius 
suggested by the gas emission is not a unique property of 
the HV Tau C disk. Indeed, this phenomenon is rather 
common among protoplanetary disks (e.g., Pietu et al. 
2005; Isella et al. 2007; Panic et al. 2009), although few 
objects have been mapped at very high resolution in both 
the continuum and line emission. Hughes et al. (2008) 
have recently demonstrated that this can be reproduced 
if the surface density profile of disks is tapered, rather 
than sharply truncated, outside of a certain radius, for 
instance. Based on the gas emission alone, the disk outer 
radius might be as large as 100 AU, twice the estimate 
from the continuum and scattered light images. Since 
our entire analysis focuses on the dust component of the 
disk, we do not consider such large disk radii in our mod- 
eling. 

4. PANCHROMATIC MODELING 

The objective of this section is to compare our broad 
set of observations of HV Tau C to predictions of a radia- 
tive transfer code in order to constrain the main prop- 
erties of the disk. Since this is the first attempt at si- 
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multaneously reproducing scattered light images, ther- 
mal emission maps and the SED of HV Tau C, our ap- 
proach consists in considering a simplified parametrized 
disk structure in an effort to search for a model that 
would represent a reasonable global fit to, rather than 
an exact representation of, all observations of this disk. 

4.1. Radiative transfer models: overall framework 

We use the MCFOST radiative transfer code (Pintc 
et al. 2006) which computes synthetic observables, such 
as SEDs and images, by propagating photon packets 
through the disk. Scattering, absorption and reemission 
by dust grains are taken into account following the Mie 
theory valid for homogeneous spherical grains. We as- 
sume that the dust grains are at the local thermal equi- 
librium with the surrounding radiation field throughout 
the disk. Synthetic temperature maps, SEDs and im- 
ages are computed simultaneously at all inclinations, al- 
lowing us to fit for that parameter as well. For each 
model, the temperature, SEDs and images computation 
are performed using typically 1, 4 and 32 million packets, 
respectively. The disk is assumed to be passive, i.e., its 
only source of heating is radiation from the central star. 

The disk geometry is described by its inner and outer 
radii, i?i n (which we fix at 0.15AU) and i? ut, as well 
as by two independent power laws. The surface density 
profile of the disk follows S(r) = Xo(r/ro) Q . Vertically, 
we assume that p(r,z) = po(r) cxp[— z 2 /(2_ff (r) 2 )], ap- 
propriate for a vertically isothermal, non-self-gravitating 
disk in hydrostatic equilibrium. Finally, we adopt a flar- 
ing law for the scale height, namely H(r) = Ho(r /r^) 13 ■ 
The dust content of the disk is characterized by a power 
law size distribution, dN(a) oc a~ p da from a m [ n to a max . 
We adopt the optical properties of the commonly-used 
"astronomical silicate" mixture from Draine (2003) with 
p = 3.7 and a mm — 0.03^m. In all our models, we adopt 
a default distance of 140 pc, consistent with the system's 
parallax estimated by Bcrtout & Genova (2006). The 
central star is described by an effective temperature of 
3800 K based on it spectral type (Woitas & Leinert 1998; 
White & Ghez 2001; Appenzeller et al. 2005), and we use 
the corresponding log g = 4.0 NextGen synthetic spec- 
trum (Baraffe et al. 1998). Although is difficult to 
constrain in HV Tau C since we cannot measure the to- 
tal bolometric luminosity of the object, we consider it as 
a free parameter since there is no simple method to fix 
it a priori. 

Stapelfeldt et al. (2003) showed that adding an enve- 
lope to the system yields a much improved fit to the 
visible scattered light images of HV Tau C, especially 
to account for the roughly spherical halo that is detected 
well above the disk midplanc. To implement this, we add 
a spherically symmetric envelope to our model, which is 
characterized by i?; n = 1 AU and i? ut = 85 AU (maxi- 
mum extent of the "rays" identified by Stapelfeldt et al. 
2003) to mimic the halo seen in the HST images. For 
this envelope, we use a total dust mass of 5.1O _8 M , 
a radial density profile p(r) oc r~°- 5 and interstellar-like 
dust grains (a max = lpm and p = 3.7). The optical 
depth through the envelope at A = 0.5pm is r « 0.5, 
i.e. the maximum acceptable considering the morphol- 
ogy of the visible images of the disk. While a spher- 
ical envelope has little physical relevance, as its free- 



fall time would be much shorter than the system age, 
we do not explore more sophisticated envelope models 
(e.g., rotation-supported). Indeed, we consider our ap- 
proach as a reasonable "placeholder" to avoid systemat- 
ically poor x 2 values for our model scattered light im- 
ages. Images of HV Tau C show that the envelope is 
only significantly detected shortward of lpm, anyway, 
and even then the region of highest signal-to-noise ratio, 
which contributes most to the total x 2 is dominated by 
scattering off the disk. We therefore do not expect our 
simplistic parametrization to induce any strong bias on 
our results. Finally, we note that the asymmetry of the 
disk seen in scattered light image produces such strong 
deviations from our simple model that attempting to ad- 
just the envelope profile is of little interest at this point. 

4.2. Parameter space exploration 
4.2.1. Strategy and modeling grid 

While we have set a number of parameters of the 
model, we still have to deal with an 8-dimensional param- 
eter space: the total dust mass in the disk (Md us t), -Rout, 
H (defined at r = 50 AU), /3, a, a max , i?*, and the in- 
clination to our line of sight (i, where i = 0° corresponds 
to a face-on disk) are all free parameters. When comput- 
ing synthetic SEDs, we also added a foreground extinc- 
tion Ay following the interstellar extinction law which 
we considered as an additional free parameter, with val- 
ues ranging from to 5 mag. This extinction represents 
attenuation by material located around or between the 
HV Tau system and us and does not include attenuation 
by the edge-on disk and spherical envelope themselves. 

For all model parameters except the inclination, we 
adopt a coarse sampling that encompasses reasonable es- 
timates either from the previous modeling by Stapelfeldt 
et al. (2003) and Andrews & Williams (2005) or our 
own estimates above (see Table 4). For instance, we 
chose a range of values for a max that extends from lpm 
to 1cm. As shown in Figure 8, a max — lpm results in 
an opacity law that is in reasonable agreement with the 
measured interstellar extinction law through the mid- 
infrared. By extension, we will refer to this model as 
"interstellar dust" . The only parameter for which we 
use a large number of possible values is the inclination. 
MCFOST produces synthetic disk observables at all in- 
clinations at once. The final images and SED are stored 
in several "inclination bins" which are equally spaced in 
cos(i), appropriate for a random orientation prior. By 
using 91 independent inclination bins in all our simula- 
tions, we obtain a 0?6 sampling in the vicinity of the 
perfectly edge-on geometry. 

Exploring the parameter space with such a coarse sam- 
pling is computationally manageable, although it comes 
at the cost of failing to find a model that fits perfectly all 
data. Our original goal was to use this grid to identify 
a small region of the parameter space that produced a 
good fit to all observations and to run a second, finer grid 
in that smaller regions. However, as we discuss below, 
our analysis revealed that no single region of the param- 
eter space could be selected that way, as different ob- 
servations point towards different parts of the parameter 
space. We discuss possible avenues for improvement for 
future modeling efforts in Section 5. With 8 free param- 
eters in our model grid, we have computed over a million 
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independent models, which required about 126,000 h of 
CPU time on a 400-processor cluster, 92% of which was 
devoted to the synthetic SED calculation. 

Because we chose a fixed value for R in , the temperature 
at the inner edge of the disk varies from 800 to 1300 K 
in the grid. Similarly, the temperature at the outer edge 
ranges from 15 to 30 K. The radial dependence of the 
dust temperature in the midplane outside of ~ 0.25 AU 
closely follows a power law, T oc r~ q , whose index is 
in the range q = 0.4-0.6 and mostly depends on the 
amount of flaring in the disk. The dust temperature at 
a given location depends primarily on L*, hence i?* and, 
to a lesser extent, on a max (via the total dust opacity) 
and other geometric parameters. Our modeling strongly 
favors R± = 3i? Q (see below), for which the range in 
maximum temperature in the disk is 1 100-1350 K, close 
to the sublimation temperature of dust grains. While 
setting R in on the basis of this sublimation temperature 
for all models may be more physically relevant (Muzerolle 
et al. 2003), it is a computationaly expensive iterative 
process, and we believe that it would not dramatically 
change our conclusions. 

We decided to model 5 independent observables of 
HV Tau C, which we consider as representative of the en- 
tire dataset available to us: its broadband SED (Table 3), 
the F814W, H and V (LGS) scattered light images, and 
our IRAM 1.3mm visibilities as a function of projected 
baselines. Our choice of wavelengths for the scattered 
light images to fit for is a trade-off between considering as 
wide a wavelength range as possible while securing data 
that has high enough spatial resolution and signal-to- 
noise per pixel. The latter criterion led us to discard our 
new M s image in this modeling. In any case, a model that 
reproduces the F814W, H and V images reasonably well 
is likely to also reproduce scattered light images at other 
wavelengths throughout the entire 0.8-5/xm range. MC- 
FOST synthetic scattered light images were computed 
as monochromatic images at the effective central wave- 
length of each filter and with a 2" total field-of-view. The 
thermal emission from both the disk and the envelope 
is neglected in computing these images, as most of the 
emission arises from the inner most regions, which is vir- 
tually unresolved as seen from the disk outer edge, where 
scattering towards the observer occurs. As discussed be- 
low, we only aim at reproducing the morphology and we 
have checked that this morphology is virtually unchanged 
if we take into account the disk emission, which other- 
wise induces a heavy computational cost. In the millime- 
ter regime, we computed 1.3 mm thermal emission maps 
with a very fine spatial sampling to avoid aliasing in the 
Fourier transform (O'.'05/pixel). 

4.2.2. Goodness-of-fit estimates 

While we are interested in finding the best possible 
model to account for the observations of HV Tau C, we 
also aim at determining the range of acceptable models 
around it. For this purpose, we adopt a Bayesian infer- 
ence method (e.g., Lay et al. 1997; Akcson et al. 2002; 
Pinte et al. 2008), in which each model is assigned a 
probability that the data are drawn from the model pa- 
rameters. In cases where the prior has a uniform prob- 
ability distribution, as is the case here either on a lin- 
ear, cosinus-like or logarithmic scale, this probability is 
P = Pq cxp(— x 2 /2), where x 2 1S the reduced chi square. 



The normalization constant, Pq, is chosen so that the 
sum of the probabilities over all models in the grid is 
unity. Once this is done for all models in our grid, we 
can derive the probability distribution for a given pa- 
rameter by marginalizing the 8-dimensional probability 
hypercube against the other 7 dimensions. 

The computation of Xsed has 12 degrees of freedom. 
To take into account the noise intrinsic to Monte Carlo 
simulations, we ran 10 independent realizations of one 
of the best-fitting model and quadratically added the re- 
sulting standard deviation to the observational uncer- 
tainties. For the number of packets we used, this Monte 
Carlo noise is about 10% in the optical/near-infrared and 
(sub)millimcter regimes but reaches about 40% in the 
mid-infrared where photons have the hardest time escap- 
ing their deeply embedded emission region. Even though 
it is likely that this overestimates the uncertainties for 
non-edge-on models (for which photons escape more eas- 
ily to the observer) , it is of little importance since these 
models have very poor x 2 values to start with, given 
that they do not produce the characteristic double- hump 
SED. In computing Xsed' we chose to fit for \n(vF v ) in- 
stead of vF v , as this better handles the case of datasets 
that are dominated by calibration uncertainties (i.e., all 
measurements beyond 100^m), which are multiplicative 
rather than additive in nature. 

For scattered light images, we perform a pixel-by-pixel 
computation, though with a few adjustments. First of 
all, we resampled the PI band image to a pixel scale of 
0'.'027 for a higher signal-to-noise without compromising 
the spatial resolution. All synthetic images were then 
convolved with the same PSFs as used for the deconvolu- 
tion of HV Tau AB (and a TINYTIM-generated F814W 
PSF 3 ). We then aligned the model and empirical images 
using cross-correlation in the Fourier domain to deter- 
mine the offsets with which the model images are best 
aligned with the observed ones. We further normalized 
all images to a peak value of unity since our goal is to 
reproduce the morphology of the images, the flux being 
fitted for in the SED. Finally, to avoid including many 
pixels where both the data and model are indistinguish- 
able from zero and that would artificially improve our x 2 > 
we selected to use square binary masks encompassing the 
area where there is significant flux from HV Tau C. These 
masks are 2", 176 and 174 on a side at F814W, H and L', 
respectively, resulting in 2017, 3713 and 2201 degrees of 
freedom. Tests conducted without these masks showed 
that this does not affect the ranking of the models from 
best to worse while ensuring that the best x 2 are not arti- 
ficially low, which would result in too relaxed constraints. 
Similar to the SED x 2 computation, we estimated Monte 
Carlo noise maps by computing 10 realization of the best- 
fitting model at all three wavelengths. For the number of 
packets we used here, relative uncertainties range from 
less than 2% at the image peak to 8-12% at the disk's 
outer edge. 

At 1.3 mm, the available data are a series of correlated 
fluxes as a function of projected baselines. To reduce the 
number of visibilities to compute for each model and to 
improve their signal-to-noise ratio, we first flagged out 
all correlated fluxes whose uncertainty was larger than 

3 http://www.stsci.edu/software/tinytim/tinytim.html 
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0.1 Jy, i.e., twive the total flux from the source itself. We 
then averaged the correlated fluxes in the (u,v) plane 
using a 50 m bin size. This yields 60 independent mea- 
surements, or 52 degrees of freedom, at 1.3 mm. The 
model images were first rotated so that the disk mid- 
plane lies at a position angle of 108?3, and convolved 
with a 0'.'3 Gaussian "seeing" . After padding the images 
with zeros to obtain the adequate total field of view, we 
computed their Fast Fourier Transform to obtain a set 
of synthetic correlated fluxes. As for the treatment of 
the scattered light described above, here we do not wish 
to fit for the total flux in the 1.3 mm map but for its 
morphology. Therefore, we convert the correlated fluxes 
into visibilities using the observed 49.5 mJy flux for the 
input data and the measured total flux in each map for 
the models. The 1.3 mm synthetic maps are essentially 
noiseless and we neglect the Monte Carlo noise in this 
regime. 

4.3. Modeling results 
4.3.1. Best-fitting models 

The x 2 values for the models that best fit each obser- 
vation of HV Tau C are listed in Table 5. The combi- 
nation of observational and numerical uncertainties used 
in computing the goodncss-of-fit of a given model is un- 
likely to follow exactly a normal distribution. Caution 
should therefore be used when interpreting individual x 2 
values as they could be systematically biased. Nonethe- 
less, the ranking of the models is likely correct and the 
model with the absolute lowest x 2 value should be close 
to what would have been the best possible fit within the 
explored parameter space in the absence of non-Gaussian 
noise. 

First of all, we note that the best reduced x 2 f° r the 
SED, H- and L'-band images fit are in the 3.7-4.9 range. 
Considering the intrinsic asymmetry of the disk and the 
coarse sampling of the parameter space, we consider that 
these observations are well reproduced. While the in- 
trinsic structure of protoplanetary disks is likely to be 
much more complex than the power law used here, the 
simplified parametrization adopted here appears reason- 
able and can therefore provide valuable insight about the 
disk properties. The fact that the best Xfsi4W ^ s onr y 
11.5, i.e., much worse than the other two scattered light 
images, results in part from its higher signal-to-noise ra- 
tio which amplifies any departure from the models. In 
addition, the relative undersampling of the WFPC2 im- 
age results in a poorer ability to register images, thereby 
globally increasing all values of Xfsi4W Lastly, we note 
that the range of Xi.3mm m our entire grid is very narrow, 
about a factor of 2 from best to worst, implying that this 
observation of HV Tau C is only mildly constraining for 
our model. This was to be expected considering that the 
disk is only marginally resolved in our PdBI data. 

Figure 9 compares the observed SED of HV Tau C to 
several "best" models. The best overall model has a fore- 
ground extinction of Ay = 1 mag, marginally lower than 
the extinction we estimated for HV Tau AB. Both the 
model that best fits the F814W image and the best over- 
all model with a max = 1/im (i.e., interstellar-like dust) 
fail badly at reproducing the photometric data, with 
Xsed °f 139.9 and 20.7, respectively. Although these 
models produce a millimeter flux that is a factor of at 



least a few lower than the observed flux, it is interest- 
ing to note that they yield a shallow millimeter spectral 
index as a result of their high optical depth. Indeed, 
all models that adequately reproduce at least one of the 
HV Tau C observations have ri. 3mm > 15 in the disk mid- 
plane. This is a clear reminder that millimeter spectral 
indices can be strongly affected by optical depth effects 
in the case of edge-on disks. Another shortcoming of an 
interstellar-like dust model, as discussed by Pontoppidan 
ct al. (2007), is that it yields a near- to mid-infrared slope 
that is far too steep compared to observations if they are 
required to produce the right flux ratio between the two 
humps of the SED. This is a result of the vanishingly 
small mid-infrared albedo of any dust model that does 
not include grains of several microns in size. Therefore, 
the SED of HV Tau C points to the presence of at least 
intermediate-size grains in the disk. 

Let us now consider the model scattered light images. 
Figure 10 compares the observed images to several key 
models. In addition to the factors listed above, the 
poorer fit to the F814W image likely arises from our 
simplistic treatment of the spherical halo, which is too 
bright at large distances above the disk midplane in the 
models. The models that best reproduce individual im- 
ages each provide a good fit to the main features of the 
HV Tau C disk, such as the distance and flux ratio be- 
tween the two nebulae. The model that offers the best 
trade-off between the three scattered light images nicely 
reproduces the wavelength dependency of d no b but has 
roughly achromatic flux ratios between the nebulae. This 
is likely a consequence of the fact that this model, which 
has a max = 1 cm, has an almost wavelength-independent 
scattering phase function, which is required to reproduce 
the intensity profiles shown in Figure 4. Interestingly, 
while the best model for each individual image corre- 
sponds to i? ut = 75 AU, the best scattered light and best 
overall model both have an outer radius of 50 AU. Again, 
the disk asymmetry prevents us from unambiguously es- 
timating the disk outer radius based on the images only. 
Finally, we note that the model that best fits the SED is 
a very poor fit to the scattered light images (x 2 ranging 
from 34.3 to 68.1). The inclination of that model is 76?3, 
giving the observer an almost grazing angle line-of-sight, 
resulting in a counternebula that is almost undetectable 
and an almost point-like appearance at L' which does 
not match the observed image at all. 

Although the 1.3 mm emission map is approximately 
reproduced by virtually all models in our grid, it is note- 
worthy that the best-fitting models have a face-on in- 
clination. This stems from the fact that the disk, if 
viewed edge-on, is unresolved along its minor axis given 
the beam of our PdBI observations. However, our data 
show a slight drop of the correlated fluxes at the longest 
baselines along that position angle. Therefore, our mod- 
eling prefers a lower inclination in order to explain this 
feature. Since we know that the HV Tau C disk is indeed 
edge-on, this may indicate that the atmospheric phase 
noise increases with baseline in a way that differs from 
our Gaussian "seeing" approximation. Alternatively, it 
could be that the intrinsic thermal emission from the 
spherical envelope is stronger than in our model, which 
only contains small dust grains that are inefficient emit- 
ters at this regime. While this could in principle induce 
a bias in our modeling, we consider this effect as negli- 
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gible considering that this dataset is the least important 
in constraining the disk model. 

4.3.2. Bayesian analysis 

Figure 11 shows the inferred probability distributions 
for each of our 8 free parameters, after marginalization 
against all 7 other parameters. Fitting for each scat- 
tered light image separately yields probability distribu- 
tions that are similar to each other and we only show 
the probability distribution for fitting all three images at 
once, for clarity purposes. The only differences between 
fits to individual images, associated with inclination and 
maximum grain size, are discussed below. The results 
from the fit to the F814W image tends to drive the results 
when grouping all three scattered light images because 
the other images are comparatively easier to reproduce, 
as shown by the lower best \u an d xh compared to the 
best Xp8i4W- 

We remind the reader that the derived probabilities 
that a given parameter takes a certain value are only 
valid within the framework of our modeling and are based 
on our approximate treatment of Monte Carlo noise in 
the simulations. Nonetheless, they represent a more reli- 
able metric in our analysis than \ 2 values for individual 
models and they therefore better highlight how each ob- 
servation informs our model of the disk. We note that 
the best individual model (shown as filled diamonds in 
Figure 11) has parameter values that correspond in most 
cases to the most likely value as determined from the 
Bayesian method, or to a value whose global probabil- 
ity is at least 20%. This demonstrates that noise in both 
data and model is well-behaved, providing support to the 
results of the Bayesian analysis. 

The scattered light images place some constraints on 
-Rout, «max, i and, to a smaller extent, M(j ust . Since 
we have used flux-normalized images, i?* remains uncon- 
strained. Similarly, a, (3 and Ho are not well constrained 
as a result of trade-off between these parameters and 
-Mdust- A scale height of 4-5 AU is preferred when fit- 
ting only the F814W but it is not conclusive and the 
other images have much flatter probability distributions 
for Hq. Our modeling of the scattered light images finds 
-Rout = 50 AU as the most likely value, although an outer 
radius of 75 AU cannot be formally excluded considering 
the disk asymmetry and total extent. The probability 
distribution for Md us t peaks at the lowest end of the 
range sampled here, and P(M dust < 10~ 4 M Q ) w 75%. It 
is worth noting that the probability distributions based 
on fitting each image individually are almost flat, but 
combining all three in a single fit yields a significant con- 
straint as it helps solving for some of the ambiguities 
between parameters. Similarly, we find that P(a max < 
100/xm) « 90% when fitting all three images at once, but 
there is an underlying trend as a function of wavelength. 
Indeed, P(a max < Ifim) drops from 45% (most likely) 
at 0.8^m, to 20% at 1.6/xm and to 7% (least likely) at 
3.8^m, respectively. This gradual and smooth evolution 
with wavelength is probably indicative of a real physi- 
cal phenomenon. Finally, the combined fit to the scat- 
tered light images suggests an inclination 4 of 82?7 + °ig. 

4 Uncertainty ranges are defined by the 34-percentile on each 
side of the most likely value. 



However, there is a marginal (2a) trend as a function 
of wavelength: the most likely inclinations are 80?8 + i i g , 

84?6+/-2?5 and 85?9+^ for the F814W, H and V im- 
ages, respectively. This trend stems from the inability of 
our model to produce a wavelength-dependent flux ratio 
between the nebulae without generating chromatic vari- 
ations of the intensity profile along the major axis. This 
may indicate that our dust model does not possess the 
adequate chromatic behavior. 

Considering now the SED, our modeling mostly pro- 
vides constraints on the same parameters as the scat- 
tered light images. The SED fits favor the presence of 
large grains and a high total dust mass, with P(a max > 
lOO^m) w 85% and P(M dust > 10~ 4 M Q ) w 85%. Both 
conditions are necessary to produce sufficient fluxes at 
the long wavelength end of the SED. While interstellar- 
like grains can produce the observed spectral index, all 
models with such a dust population fall short in the 1- 
3 mm range by at least an order of magnitude. While 
SED fitting is generally not sensitive to the outer radius 
of a disk, here we find that a small outer radius is strongly 
preferred: P(R out < 50AU)/P(R out > 75 AU) w 4. 
This is because the relative height of the two humps in 
the SED and the depth of the trough that separates them 
is directly influenced by the disk geometry, via the total 
optical depths along our line of sight to the star on one 
hand and from the star to the upper layers of the outer 
disk on the other hand. Finally, the modeling of the 

SED yields a best fit inclination of 75?0 +4 i 5 . The fact 

J -3.9 

that the SED provides a weaker constraint on the disk 
inclination than the scattered light images is typical of all 
disk analyses. Nonetheless, the particular viewing angle 
of edge-on disk allows for a reasonably narrow range of 
inclination to be defined, and it is reassuring to note that 
the difference between the inclinations derived from the 
scattered light images and the SED is at the insignificant 
1.6a level. 

The fit to the 1.3 mm yields much weaker constraints 
on the model parameters: all associated probability dis- 
tributions are consistent with being flat except for i? ut, 
a and, to a smaller extent, i. As discussed above, this is 
a result of the fact that the disk is only marginally re- 
solved in our PdBI data. To best reproduce the observed 
visibilities, the model favors a face-on geometry, a large 
outer radius and a flat surface density that places a lot of 
mass, hence of millimeter emission, at large radii. As far 
as inclination is concerned, the inferred most likely value 
is +50 °, illustrating the weak preference for a face-on 
geometry. The constraints on both R out and a are much 
stronger with a ratio of the most- to least-likely value of 
at least 4:1. 

The strength of our modeling approach is to be able 
to fold all independent datasets in a single coherent fit, 
which both yields much sharper constraints and a view 
of the intrinsic contradictions between observations. For 
some parameters, such as (3, Hq and R+, the combined 
fit provides better constraints than any individual obser- 
vation as a consequence of the fact that each observa- 
tion is associated to a different set of ambiguities. We 
thus find that the disk is significantly flared (J3 > 1.15), 
that its scale height at 50 AU is 5AU at most, and that 
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R± ~ 3i?0, which implies ~ 1.7X0. Assuming that 
the accretion luminosity is negligible, this is consistent 
with a lMyr 0.7-1 M central star (Baraffe et al. 1998), 
sightly larger than the mass estimated for HV Tau A 
(White & Ghez 2001). Although uncertainties are large, 
this is in rough agreement with the dynamical mass that 
we have estimated from the rotating CO disk. 

In the case of i? ut and a, the constraints set by the 
millimeter mapping are superseded by those from the 
scattered light images and the SED because the former 
is comparatively too easy to reproduce in our grid. While 
a 1/r surface density profile is preferred in our model, all 
other values tested in our grid have non-negligible prob- 
abilities. Similarly, our overall fit favors i? ut = 50 AU 
only by a factor of 2 over the larger radius. Excluding 
the 1.3 mm map from our analysis increases the prefer- 
ence for the smaller outer radius to a factor of 4:1, how- 
ever. Future higher resolution millimeter observations 
are needed to better constrain both of these parameters. 
Also, the final probability distribution for the disk in- 
clination is largely driven by the scattered light images. 

, 2° c 

Our global fit suggests an inclination of 82?1 

The most interesting results of our analysis concern 
flmax and Mdust, for which scattered light images and 
SED lead to contradictory predictions (see Figure 11. 
The latter favors large grains and a high dust mass 
whereas the former suggests small grains and a low dust 
mass. Our global fit, which aims at finding the best pos- 
sible trade-off between all constraints points to a very 
high dust mass (Md us t > 10~ 3 M Q ) and an intermediate 
maximum grain size (a max w 10/im). Since both of these 
values are essentially rejected (P < 10%) by at least 
one of the observations, this "best" model should not be 
considered as a good fit. Rather, this apparent contra- 
diction between the various observations, which probably 
mirrors the trend in preferred a max as a function of wave- 
length for the scattered light images, indicates that our 
model is too simplistic in at least some aspects. Possible 
explanations are explored in Section 5. 

4.4. Comparison to previous modeling 

There are three published models of the HV Tau C 
disk. In their discovery study, Monin & Bouvier (2000) 
derived a disk inclination of w 84° based on a single scat- 
tering model. That inclination has later been confirmed 
by subsequent models of the disk, including the present 
work. The second, most sophisticated to date, model of 
the disk was put forth by Stapelfeldt et al. (2003) who 
analyzed scattered light images at 0.8 and 2.2/zm. More 
recently, Andrews & Williams (2005) derived a total 
disk mass for HV Tau C based on its sub-millimeter flux 
although they did not conduct a complete SED fit for 
this source. Here we compare our model results to these 
previous efforts. 

The model constructed by Stapelfeldt et al. (2003) was 
based on the same F814W image we have used here, and 
a lower quality 2.2^m adaptive optics image that was 
characterized by a poorer resolution (0V13 compared to 
0'.'08 for our VLT/NAOS image) and a smaller achieved 
Strchl ratio. Their conclusions are similar to ours in 
several ways. First of all, they found that it is easier 
to fit for the K band image, as demonstrated by the 
lower achieved \ 2 value. Second, fitting for the envelope- 



free K band image alone points towards an inclination 
that is closer to edge-on than if the F814W image is 
included in the fit. Conversely, when fitting both images 
simultaneously, they also note that they do not succeed 
in reproducing the wavelength dependence of the flux 
ratio between the nebulae. 

If we temporarily focus on the F814W image alone, 
which is the most constraining dataset in their fit, our 
Bayesian analysis yields a max = 1/im, Md us t = 10 _5 M Q 
and H = 4AU and i = 80? 8. These values are con- 
sistent with those of Stapelfeldt et al. (2003) within our 
uncertainties. The derived inclination, somewhat fur- 
ther away from edge-on, may be explained by our differ- 
ent treatment of the spherical halo, which plays a non- 
negligible role in shaping up the morphology of HV Tau C 
at 0.8/zm. Finally, we note that our simultaneous fit 
to the F814W, H and V images favors a max = 1/zm. 
This model has a ratio of opacity of Ko.s^m/^^jum ~ 3.3 
and scattering asymmetry parameters of go.s^m = 0.6 
and <72.2/jm = 0.55. These properties are very similar 
to those derived by Stapelfeldt et al. (2003), namely 

K0.8ftm/ K2.2nm ~ 3.5 and go.&fim = 0.65. 

Overall, we find that our fit to scattered light images is 
in good agreement with the previous modeling effort for 
this source, despite differences in strategy and parame- 
ter space covered. However, the models that best repro- 
duce the F814W image produce millimeter fluxes that 
are one to two orders of magnitude too low in the mil- 
limeter regime, in apparent contradiction with the con- 
clusion of Stapelfeldt et al. (2003), who found agreement 
within a factor 2. The explanation for this discrepancy 
is that these authors assumed a ratio of dust opacity of 
^o. 8^111/^1. 3mm = 6.10 3 . However, with our assumed dust 
composition, this ratio is about 7.10 4 for a max = 1/im, 
accounting for most of the difference in predicted mil- 
limeter flux between the two models. 

Andrews & Williams (2005) presented submillimeter 
observations of HV Tau C, which we have included in our 
fit. Because the source showed an unusual spectral index, 
they did not attempt to fit a full-fledged model but rather 
used an empirical recipe to convert the measured 850^m 
flux into a disk mass. They derived Md us t = 2.10 -5 M , 
almost two orders of magnitude lower than the value we 
have derived here based on cither the SED or global fit. 
The overall shape of the SED for HV Tau C suggests 
that the 850/zm flux is underestimated, possibly a fac- 
tor of 3 or so, accounting for part of this discrepancy. 
In addition, the empirical law that Andrews & Williams 
(2005) used was mostly derived from sources which are 
not edge-on and for which optical depth effects are negli- 
gible. However, the model defined by the combination of 
each of the most likely parameter value has r 1.3mm ^ 100 
in the disk midplane. We believe that this factor is the 
most important in explaining the difference in the disk 
masses inferred here and by Andrews & Williams (2005). 

4.5. Shortcomings of the model 

While our model has been successful at reproducing 
the SED of HV Tau C, as well the morphology and sev- 
eral key chromatic behavior of its scattered light images, 
it falls short in several aspects which are worth explor- 
ing. As expected, a significant shortcoming of our model 
is its built-in axisymmetric assumption, which prevents 
us from finding perfect matches to any of the scattered 
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light images. Since the nature of the asymmetry can 
only be speculated upon, no clear path to resolution can 
be provided until further monitoring clarifies its origin. 
This peculiarity is however not sufficient to account for 
the poor Xtot 01 our overall best fitting model. The most 
glaring limitation of our model is the fact that it can- 
not account simultaneously for the SED on one hand 
and the scattered light images on the other. This fact, 
which most likely is independent of the axisymmetry of 
the model, implies that more complexity has to be in- 
cluded in the model. 

The physical properties of dust grains are subject to 
improvement. Based on our modeling experience, the 
necessary requirements for an improved dust model in 
HV Tau C would be an interstellar-like opacity law and 
achromatic scattering phase function from 0.8 up 5/zm 
and an optical-to-millimcter opacity ratio < 10 4 . For 
one, one could devise a dust model that contains a frac- 
tion of C-rich grains to better conform to interstellar 
abundances. Departures from a simple power law size 
distribution are also likely. Such a phenomenon has 
been observed in the interstellar medium (Kim et al. 
1994; Weingartncr & Draine 2001) and multi-modal dis- 
tributions are predicted by dynamical models of grain 
growth, migration and fragmentation in protoplanetary 
disks (e.g., Dullcmond & Dominik 2005; Laibe et al. 
2008; Zsom & Dullemond 2008). Another direction to 
explore is the possibility of non-spherical grains. Re- 
cently, Kimura et al. (2003) successfully reproduced the 
wavelength-independent anisotropic scattering observed 
in cometary dust by considering aggregates of small par- 
ticles. Since the dust in HV Tau C also has achromatic 
scattering properties, aggregate grains are a candidate 
that deserves consideration in the future. In this context, 
it is interesting to note that fractal grains are expected 
to have a much flatter optical-to-millimcter opacity law 
(Wright 1987), akin to our models with CL max > lOO^rn, 
which are preferred in our analysis. 

Since it appears difficult to simultaneously match all of 
these constraints with a single dust model, it will likely 
be necessary to drop the assumption that the dust prop- 
erties are uniform throughout the disk. Since each type 
of observation probes a different region of the disk, de- 
coupling the dust properties in various regions of the disk 
would likely provide sufficient leeway. In this picture, 
the disk would be quite massive and would contain large 
grains, up to millimeter-sized particles, in the midplanc 
in order to explain the observed SED. However, grains 
in the surface layers would not exceed a few microns in 
size to explain the scattered light images. To limit the 
number of free parameters, we have only considered fully 
mixed models here but we suggest that this hypothesis 
be further studied in future modeling of HV Tau C. 

Future modeling efforts on HV Tau C should consider 
these various possibilities in improving on the model pre- 
sented here. A robust strategy to move forward could 
consist in conducting detailed modeling of a particular 
type of observations, such as that conducted by Wat- 
son & Stapelfeldt (2004) on the scattered light images 
of HH 30 for instance. These specific analyses would al- 
low one to extract the maximum amount of information 
from each observation and to feed it into a refined global 
model. 



5. DISCUSSION 
5.1. Global properties of the HV Tau C disk 

Our objective in this work was to conduct the first 
"global" fit to the structure and dust content of the 
HV Tau C disk, one of the few objects in which such 
an effort is conducted in a self-coherent way. Although 
some uncertainty remains, our analysis has led to some 
robust conclusions about the HV Tau C disk. For in- 
stance, the total mass of the disk must be at the high 
end of the range probed here. Even in the presence of 
large grains, which are efficient emitters in the millime- 
ter regime, the total dust mass has to be on the order 
of Mdi s k ~ 10~ 3 M Q . Indeed, such a high dust mass 
is necessary to account for the observed fluxes and rel- 
atively high optical depths up to 3.3 mm. Assuming a 
canonical 100:1 gas-to-dust mass ratio, this implies that 
the disk is quite massive compared to the central star 
(M disk /M* > 0.1). This is at the high end of the dis- 
tribution observed for disks surrounding TTauri stars 
(Andrews & Williams 2007). It is plausible that the disk 
is only marginally stable, in which case the presence of 
spiral density waves could be responsible for the observed 
asymmetry and/or variability of the system. 

One respect in which the HV Tau C disk is special 
is its compact radius, a probable consequence of tidal 
forces induced by HV Tau AB. With R out = 50 AU 
and a = —1, we derive a total surface density of 2800 
g cm~ 2 1 AU away from the star, comparable to the sur- 
face density inferred for the early Solar Nebula as well 
as for extra-solar planetary systems, albeit with a pos- 
sibly flatter surface density profile (Hayashi 1981; Kuch- 
ner 2004). Indeed, the surface density at large radii in 
the disk, almost 60g.cm~ 2 at 50 AU, appears substan- 
tially higher than those derived for other protoplanetary 
disks (Dutrey et al. 1996; Kitamura et al. 2002). The 
disk around HV Tau C is a clear example of the fact 
that binary systems can host circumstellar disks massive 
enough to form planets for timescales of several million 
years. Indeed, HV Tau may well be a prototype for the 
earliest evolutionary stages of field stars that are found 
to host both an extrasolar planet and at least one stellar 
companion (e.g., Eggcnberger et al. 2007; Mugraucr & 
Neuhauser 2009). 

As a test of the self-consistency of our model, it is 
interesting to compare the disk scale height we have 
derived, Ho < 5AU, to the simple assumption of ver- 
tical hydrostatic equilibrium that is embedded in our 
parametric disk structure. Assuming molecular hy- 
drogen, the scale height can be written as iJQ ydro = 
3A(T/20K) 1 / 2 (R/50AU) 3 / 2 (M i ,/M & )~ 1 / 2 AU. Our best 
global model has a midplane temperature of 19 K at 
50 AU, so that our upper limit on Hq implies M* > 
0.44M Q . This is not a stringent constraint on the mass of 
the central star, but we note that this lower limit does not 
violate the mass inferred from the spectral type and lu- 
minosity of HV Tau C nor the rough kinematic estimate 
we derived from the rotating CO disk. Both estimates 
are in the 0.7-1 M Q range, for which the derived scale 
height would be 3.3-4 AU at 50 AU, which is consistent 
with our modeling. We therefore conclude that the outer 
disk can reasonably be described by a disk in hydrostatic 
equilibrium. In addition, we suggest that HV Tau C may 
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be the most massive component in the triple system. 

Finally, while large uncertainties remain regarding the 
dust properties, it is clear that dust grains much larger 
than l\im are found in at least some parts of the disk. 
This could be evidence for grain growth in the disk, al- 
though it is plausible that growth to a few microns may 
have predated the formation of the disk itself (McCabe 
et al. 2003). Indeed, analyses of extinction laws have 
already pointed to the presence of small quantities of 
grains a few microns in size in the interstellar medium 
and even more so in molecular clouds (Kim et al. 1994; 
Weingartner & Draine 2001; Indebetouw et al. 2005; Fla- 
herty et al. 2007). On the other hand, growth to mm- 
sized grains, which seems supported by the analysis of 
the SED of HV Tau C, has most likely occurred in the 
disk itself. If the disk is indeed vertically stratified and 
no such grains are found in the disk surface, it remains 
unclear whether this is because of a very efficient settling 
process or a mere consequence of the difficulty to grow 
such large particles in the lower density regions of the 
disk. 

5.2. Comparison to other circumstellar disks 

As outlined in the introduction, evidence for grain 
growth and dust stratification in protoplanetary disks 
has been mounting over the last two decades. It is there- 
fore not surprising that our modeling of the HV Tau C 
disk calls for both processes. However, one must keep 
in mind that a given observation only probes a limited 
region of the disk. For instance, mid-infrared emission 
is dominated by emission from the inner few AUs at 
most, providing no information on dust stratification in 
the outer disk. In addition, there is an intrinsic observa- 
tional bias in all studies whereby only a small range of 
grain sizes, those with the largest effective cross-section, 
actually contribute to the observed fluxes. Therefore 
longer wavelength observations tend to systematically 
call for larger grain sizes, even if one limits the analy- 
sis to scattered light images (e.g., McCabe et al. 2003). 
Similarly, the mid-infrared silicate features arise from mi- 
cronic and submicronic grains, providing no constraint 
on much larger particles. To draw a complete picture of 
dust properties throughout a given disk, it is therefore 
necessary to gather wide and homogeneous datasets, as 
we have done here. 

As we have already mentioned, few objects have been 
studied in as much detail as our present analysis of 
HV Tau C. A clear dichotomy in terms of grain size 
was also found in the case of 'iRAS 04302+2247 (Wolf 
et al. 2003) but the interstellar-like dust grains favored 
by the scattered light images were located in the rela- 
tively massive envelope surrounding this embedded ob- 
ject and not in the disk itself. Vertical stratification was 
demonstrated in the massive GG Tau ring (Duchene et 
al. 2004), although the available evidence only applies 
to grain smaller than 10/im. The HKTauB edge-on 
disk also appears to have a vertically stratified structure 
(Duchene et al. 2003), although a self-consistent mod- 
eling of its SED, thermal emission maps and scattered 
light images is still required to make direct comparison 
with HV Tau C. On the other hand, Glauser et al. (2008) 
were able to reproduce a wide range of observations of 
IRAS 04158+2805 without the need to include vertical 
stratification. A similar conclusion was recently reached 



by Sauter et al. (2009) for CB 26. Overall, these detailed 
analyses confirm the general trends outlined above that 
both grain growth and stratification are common, but not 
ubiquitous, in protoplanetary disks. They further offer 
the unique opportunity to test the detailed physics of 
dust evolution, such as the amount of dust stratification 
or its radial dependency (e.g., Pinte et al. 2007). 

An interesting head-to-head comparison can be drawn 
between HV Tau C and IMLup, for which Pinte et al. 
(2008) recently conducted a similar panchromatic mod- 
cling analysis, including scattered light images, millime- 
ter emission maps and the system's SED. The central 
stars have the same spectral type, and similarly high disk 
masses were inferred, making these two systems very sim- 
ilar. Contrary to the conclusion reached here however, 
Pinte et al. (2008) were able to place stringent constraints 
on all disk parameters and one may wonder what the rea- 
sons for this difference are. For one, the IMLup disk is 
not seen edge-on, which allows the authors to fix some of 
the model parameters (e.g., i?*, dust composition, M^j in 
an unambiguous way. The smaller dimensionality of the 
parameter space to explore allowed Pinte et al. (2008) to 
sample it more finely. Another key feature of the IM Lup 
disk is that it is much less compact than HV Tau C 
and is well resolved at millimeter wavelengths, providing 
more stringent constraints on the model. Last but not 
least, a qualitative analysis of the SED of IMLup read- 
ily demonstrated the need for a stratified structure in the 
disk, which introduces an extra degree of freedom that we 
could not afford in our analysis of HV Tau C. Nonethe- 
less, our study has shown that significant constraints on 
the disk properties can be obtained in an edge-on geom- 
etry. No global modeling of HH30 or HKTauB, other 
prototypical edge-on disks, are yet available, but con- 
ducting such efforts should yield some interesting results 
and are needed to move forward in our understanding of 
protoplanetary disks. 

6. CONCLUSION 

We have obtained new l-5^im adaptive optics images 
of the HV Tau triple system at the VLT and Keck obser- 
vatories, including the first 4.8/xm scattered light image 
of the edge-on HV Tau C disk. All of our images reveal a 
steady lateral asymmetry in the disk that indicates a de- 
parture from pure axisymmetry in the disk. This could 
be related to the known variability of HV Tau C. We ex- 
tact precise relative astrometry from our new images and 
find that HV Tau AB-C constitutes a common proper 
motion pair. We also resolved the tight binary system 
HV Tau AB and found a surprisingly slow orbital mo- 
tion compared to its projected separation, probably due 
to orbit eccentricity and/or a large deprojected physical 
separation. We have also obtained the first spatially re- 
solved 1.3 mm continuum and 12 CO 2-1 millimeter maps 
of the HV Tau C disk with IRAM's PdBI interferome- 
ter. The continuum emission is resolved along the same 
position angle as the scattered light images, as expected 
from dust thermal emission. The CO map shows evi- 
dence of Keplerian rotation about an 0.5-lM Q central 
star. We also completed the SED of that component 
with new Spitzer/MIPS 24, 70 and 160/im observations 
as well as observations at 350/im, 2.8 and 3.3 mm at the 
CSO, PdBI and CARMA facilities. The 1-3 mm spectral 
index of HV Tau C is relatively flat (a mm = 2.5 ± 0.2), 
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indicative of the presence of millimeter-sized grains in 
the disk and/or high optical depth even at millimeter 
wavelengths. 

To interpret these data along with previously published 
HST images of HV Tau C, we have computed a grid ra- 
diative transfer models that produced synthetic disk ob- 
servations. While we can reproduce most observational 
properties of HV Tau C, a Bayesian analysis of our model 
grid reveals that the SED and scattered light images of 
HV Tau C provide mutually exclusive constraints. As 
found previously, the scattered light images are best fit 
with a relatively small total mass and interstellar-like 
dust properties. Fitting the SED, however, requires al- 
most two orders of magnitude more mass as well as a 
grain size distribution that extends to millimeter sizes. 
The mismatch between the two sets of constraints re- 
veals an intrinsic shortcoming of our parametrized model, 
which assumes perfectly mixed dust. Indeed, the small 
maximum dust grain size inferred from fitting the scat- 
tered light images is impossible to reconcile with the 
long- wavelength end of the SED. In turn, this suggests 
that both grain growth and vertical stratification are 
present in the HV Tau C disk. While these phenom- 
ena have been suggested for other protoplanctary disks 
in the past, HV Tau C is one of a handful of objects 
that have been studied in sufficient details to eventually 
provide a complete picture of these processes. 

Future observations of HV Tau C can help refine the 
model proposed here. For instance, high-resolution (0'.'2 
or better) millimeter mapping of the disk would better re- 
solve it and provide more stringent constraints. In partic- 
ular, obtaining observations at even longer wavelengths, 
7 mm or 1.3 cm, would probe a regime in which the opti- 
cal depth through the disk is much smaller and enable a 
more direct interpretation in terms of disk properties. On 
the longer run, mapping of the disk with ALMA in the 
submillimcter regime and in scattered light at 8~10/um 
with the next generation of large ground-based telescopes 
will help bridge the gap between the scattering and ther- 
mal emission regimes at a roughly constant spatial res- 
olution. In particular, very high resolution maps with 
ALMA may resolve this and other disks along the ver- 
tical axis, in order to probe their vertical stratification 
and, ultimately, to outline an evolutionary sequence. 
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Fig. 1. — Scattered light images of HV Tau: the left most panel shows the F814W image obtained by Stapelfeldt et al. (2003) while the 
other three panels show our new images obtained at H and K s with VLT/NACO (the K B image is the coronagraphic one), and at L' with 
Keck/NIRC2 and the LGS adaptive optics system. The F814W image is shown on the logarithmic stretch to better highlight low-level 
features. The arrows point to the "rays" identified by Stapelfeldt et al. (2003). All other images are shown on a square root stretch. In all 
cases, the field-of-view is 6" across. As discussed in the text, the L' image has been smoothed and resampled to improve the signal-to-noise 
ratio without degrading its spatial resolution. 

Fig. 2. — Orbital motion of HV Tau B with respect to HV Tau A, whose position is indicated by the diamond. The 1994 and 1996 points 
reflect astrometric measurements by Simon et al. (1996) and Monin & Bouvicr (2000), whereas the other two points are our 2002 and 2004 
measurements (see text). The insets represent our H and L' images after a light Lucy deconvolution (25 and 50 iterations, respectively). 
While deconvolution better highlights the fact that the tight binary is resolved in our data, relative astrometry and photometry information 
for all our images was extracted by PSF-fitting of the original images. 

Fig. 3. — Contour plots of HV Tau C after rotating all frames to a common orientation where the disk midplane is horizontal. The K s 
image is the VLT/NACO coronagraphic one and the V image is the Kcck/NIRC2 LGS image. The contours lie from 90% of the peak and 
by decreasing factors of 2 from there. The L' and M B images have been smoothed to improve their sensitivity (see Section 2.1.2 for more 
details). The dashed lines are guidelines indicating the location of the peak of each nebula in the F814W image. 

Fig. 4. — Lateral intensity profiles along the brightest (NE) nebula of HV Tau C. Profiles are drawn for the F814W (solid blue), H 
(dotted green), coronagraphic K a (dashed cyan), LGS L' (dot-dashed orange) and M s (triplc-dot-dashed red) images of the disk. All 
profiles arc drawn normalizing the peak intensity in each image to unity. The thin curves represent the azimuthally-avcragcd profile of the 
PSF corresponding to each observation, using the same linestyles and color coding. The nebula is well resolved and its intensity profile is 
remarkably constant across the wavelength range probed here. 

Fig. 5. — PdBI 1.3mm correlated fluxes for HV Tau C as a function of projected baseline length for all baselines that are within 22? 5 of 
the disk major axis (squares) and within 22?5 of its minor axis (diamonds). The major axis is assumed to be 108?3 as determined from the 
scattered light images. The dashed and dotted curves represent the major- and minor-axis of the Gaussian fit to all visibilities, respectively, 
assuming a Gaussian 0''3 atmospheric "seeing". 

Fig. 6. — Top: Contour plots of the blue (dotted contours, 0.8-6 km. s -1 velocity range) and red (dashed contours, 6-9.5 km. s -1 velocity 
range) parts of the 12 CO (2-1) line, superimposed on the contours of the adjacent 1.3 mm continuum (solid contours). Contours are drawn at 
22.5, 45 and 90% of the peak intensity in each map. The reconstructed beam is shown for comparison. The origin of the relative coordinates 
is at 04h38m35.51s, +26°10'41''5 (J2000), the nominal position of HV Tau C. The inset represents the integrated line profile. Bottom: 
Position- velocity diagram along the disk major axis (indicated by the two solid segments in the top panel). The three-dimensional datacubc 
has been smoothed using a 1 km.s -1 X 1" running boxcar function. The contours in the position-velocity diagram are at 15, 30, 60 and 90% 
of the peak intensity. The dashed curve represent the theoretical Keplerian rotation curve for an 0.7 Mq star with v sys = 5.75 km. s -1 . 
This is not meant as a fit but a reference to guide the eye. The dotted lines indicate the velocity range in which 13 CO emission from the 
Taurus molecular cloud was observed by Mizuno et al. (1995) at the location of HV Tau. 

Fig. 7. — SED of HV Tau AB (asterisks) and HV Tau C (diamonds). Filled diamonds and upper limits indicate the photometry datasct 
considered in our model fitting (sec Table 3) whereas empty diamonds represent photometry at other epochs, illustrating in particular the 
optical and near-infrared variability of HV Tau C. 

Fig. 8. — Extinction law for our dust models assuming a max = lfim (thick solid), a max = 100/im (thick dotted), and a max = 1cm (thick 
dashed). Diamonds represent the measured extinction law from Cardelli et al. (1989) and Indebetouw et al. (2005). The thin solid curve 
represents a model with a m ax = 0.5/im(not used in our modeling) that best reproduces the measured interstellar extinction law up to 2fim, 
but falls short of it in the 5-8/im regime. The a m ax = l/^m is better in that regime and we use it as a proxy for "interstellar dust" in our 
model. 

Fig. 9. — Comparison of the SED of HV Tau C (from Table 3, see also Figure 7) with various models: the best fit to all available datascts 
(black solid curve), the best fit to the SED alone (red dashed curve), the best fit to the F814W image (green dot-dashed curve), and the 
best overall fit with a m ax = 1/im, i.e. interstellar-like dust properties (blue dotted curve). The thin errorbars indicated at the top of the 
plot represent our Monte Carlo noise estimate. This noise is apparent in the jumps seen at 3/^m for the best overall model and at 6/^m for 
the best F814W fit. 

Fig. 10. — Contour plots of HV Tau C at 0.8/im (top row), 1.6/im (middle row) and 3.8/im (bottom row). From left to right, the 
columns represent: the actual data, the best model in our grid at each wavelength (a), the best model to all three scattered light images 
simultaneously (6), the best overall model (c), and the best fit to the SED alone (d). All model images have been convolved with the 
appropriate PSF. In each plot, the contours lie at 80, 40, 20, 10 and 5% of the peak. 

Fig. 11. — Bayesian inference probability distributions for the free parameters in our model, based on fitting the three scattered light 
images at once (dashed blue histograms), the SED (dot-dashed red), the spatially resolved 1.3 mm correlated fluxes (dotted green), and all 
observations at once (solid black). The probability distribution for the inclination based on the 1.3 mm is almost flat and very close to 0. 
For reference, the filled diamonds indicate the value of the parameters for the model with the absolute lowest Xtot * n our grid. Notice the 
different vertical scale in the two rows of plots. 
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TABLE 1 

AsTR0METRIC a AND PHOTOMETRIC PROPERTIES OF HV TAU 



A 


Epoch 


PAB-C 


PAab-C 


PA-B 


PA A -B 


Am A -B 


w 




(") 


(°) 


(") 


(°) 


(mag) 


0.80 


2000.19 


4.04±0.01 


43.4±0.3 








1.27 


2002.90 


4.08±0.03 


43.5±1.0 








1.66 


2002.90 


4.06±0.03 


43.6±1.0 


0.054±0.003 


309±3 


0.88±0.02 


2.18 


2002.90 


4.03±0.03 


43.6±1.0 


0.062±0.004 


317±2 


0.88±0.05 


3.78 


2002.95 


4.03±0.01 


44.9±0.1 


0.063±0.003 


312±2 


0.53±0.04 


3.78 


2004.84 


4.02±0.01 


44.4±0.1 


0.058±0.002 


317.5±1.5 


0.59±0.05 


4.67 


2002.95 


4.07±0.02 


44.7±0.2 









a Astrometric uncertainties include both measurement (including centroiding) and absolute calibration uncertainties. 



TABLE 2 











Observed p 


roperties of 


the HV Tau C disk 


A 


PA disk 


<^neb 


-F^pcak 


FRint 




Note 


w 


(°) 


(") 






(") 




0.80 


108.2±1.2 


0.335±0.005 


4.0±0.2 


2.8±0.1 


1.55±0.03 




1.27 


107.5±1.1 


0.281±0.002 


2.5±0.2 


2.4±0.1 


1.04±0.07 




1.66 


108.4±0.9 


0.284±0.002 


2.3±0.2 


2.2±0.1 


0.96±0.05 




2.18 


108.4±0.7 


0.275±0.003 


1.8±0.1 


1.80±0.05 


1.21±0.07 


direct image 


2.18 


107.7±0.8 


0.267±0.003 


2.0±0.1 


1.82±0.05 


1.24±0.07 


coronagraphic image 


3.78 


108.8±1.3 


0.237±0.004 


1.9±0.2 


1.8±0.1 


0.66±0.03 


NGS image 


3.78 


109.1±0.9 


0.236±0.003 


2.1±0.2 


1.5±0.1 


0.87±0.05 


LGS image 


4.67 


108.2±1.7 


0.231±0.006 


1.8±0.3 


2.3±0.2 


0.45±0.05 





TABLE 3 
Adopted SED for HV Tau C 



A 


F v 


a(F v ) 


Date 


Ref. 


A 


Fv 


<t(F„) 


Date 


Ref. 


(H 


(mjy) 


(mjy) 






w 


(mjy) 


(mjy) 






0.545 


0.446 


0.097 


1991 Jan-Feb 


1 


11.8 


21.8 


3.3 


2002 Nov 13 


3 


0.638 


1.44 


0.31 


1991 Jan-Feb 


1 


24. 


289. 


19. 


2005 Feb 28 


5 


0.797 


2.59 


0.42 


1991 Jan-Feb 


1 


70. 


880 


63 


2005 Feb 28 


5 


1.22 


7.38 


0.97 


1996 Aug 26 


2 


160. 


<750. 




2005 Feb 28 


5 


1.65 


9.92 


0.97 


1996 Aug 26 


2 


350. 


370. 


30. 


2008 Jan 28 


6 


2.19 


8.66 


0.94 


2001 Dec 1 


3 


450. 


<519. 




2000 Feb 2 


7 


3.45 


8.44 


0.37 


2001 Dec 1 


3 


850. 


47. 


6. 


2000 Feb 2 


7 


3.6 


8.25 


1.25 


2004 Mar 7 


4 


1300. 


40. 


6. 


1988 Apr 29-30 


8 


4.5 


9.09 


0.50 


2004 Mar 7 


4 


2763 


7.1 


0.9 


2005 Feb 26 


6 


5.8 


9.39 


0.51 


2004 Mar 7 


4 


3252 


3.8 


1.2 


2008 May 29 


6 


8. 


11.4 


0.5 


2004 Mar 7 


4 













References. — 1 - Magazzu & Martin (1994); 2 - Woitas & Leinert (1998); 3 - McCabe et al. (2006); 4 - Hartmann et al. (2005); 5 - 
Spitzer Taurus Legacy Survey (Rebull et al., in prep.); 6 - this work; 7 - Andrews & Williams (2005) ; 8 - Beckwith et al. (1990) 

TABLE 4 
Model parameters 



Parameter Min. value Max. value -/V samp i Sampling 



Radiative transfer parameters 



M dust (M ) 


10" 5 


10" 3 


5 


logarithmic 


Rin (AU) 




0.15 




fixed 


i?out (AU) 


50 


75 


2 


linear 


H a (AU) 


3. 


7. 


5 


linear 





1.05 


1.25 


3 


linear 


a. 


-1.5 


0. 


4 


linear 


a min ((im) 




0.03 




fixed 


<Jmax (M m ) 


1. 


10 4 


5 


logarithmic 


P 




3.7 




fixed 


R* (Re) 


1.5 


3. 


4 


linear 


T, (K) 




3800 




fixed 




0. 


90. 


91 


linear in cos i 


Post-processing parameters 


D (pc) 




140 




fixed 


A v (mag) 





5. 


21 


linear 


^Hq is defined at 


a radius t-q = 50 AU from the central star. 
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Duchene et al. 



TABLE 5 

GOODNESS-OF-FIT ESTIMATES FOR SELECTED MODELS 





Xf814W 


xi 


xb 


Xo. 8-3. Sum 


XsED 


Y 2 

A 1.3mm 


Xtot 


Best F814W model 


11.47 


9.24 


17.69 


38.40 


139.89 


15.30 


193.59 


Best H model 


22.32 


3.74 


11.22 


37.28 


204.58 


16.62 


258.47 


Best L' model 


36.64 


5.97 


3.84 


46.45 


126.97 


14.53 


187.96 


Best images model 


16.39 


5.17 


5.80 


27.36 


174.95 


15.87 


218.17 


Best SED model 


68.10 


34.30 


42.45 


144.85 


4.94 


16.10 


165.89 


Best 1.3 mm model 


122.60 


58.80 


58.46 


239.86 


297.23 


9.19 


546.28 


Best overall model 


15.47 


7.29 


10.00 


32.76 


11.94 


16.35 


61.04 



